Show stations on the mesh¶

Read in stations info¶

In [1]:
import schimpy

from schimpy import station
dfs = station.read_station_in('../tests/data/m1_hello_schism/station.in')
dfs = dfs.reset_index()

Read in mesh info¶

In [2]:
from schimpy import schism_mesh

grid = schism_mesh.read_mesh('../tests/data/m1_hello_schism/hgrid.gr3')
In [6]:
import holoviews as hv
hv.extension('bokeh')
from holoviews import opts, dim, streams
from holoviews.operation import datashader
import warnings
warnings.filterwarnings('ignore')
import panel as pn
pn.extension()

Show points on top of mesh¶

In [7]:
trimesh = hv.TriMesh((grid.elems, grid.nodes))
# rasterize to view faster, zoom in to clarify features
img = datashader.rasterize(trimesh.edgepaths).opts(opts.Image(cmap=['darkblue'])).opts(width=800, height=400)
# spread image pixels to see mesh in a more bold style
elems_only = datashader.spread(img)

nodes_only = datashader.dynspread(datashader.rasterize(trimesh.nodes).opts(opts.Image(cmap=['blue'])), shape='circle', max_px=6)

full_mesh = elems_only*nodes_only

elems_only.opts(alpha=0.2)*hv.Points(dfs, kdims=['x','y'], 
                                     vdims=['z','id','subloc','name']).opts(color='red', 
                                                                        size=10,
                                                                        tools=['hover'])
Out[7]:
In [8]:
trimesh = hv.TriMesh((grid.elems, grid.nodes))
# rasterize to view faster, zoom in to clarify features
img = datashader.rasterize(trimesh.edgepaths).opts(opts.Image(cmap=['darkblue'])).opts(width=800, height=400)
# spread image pixels to see mesh in a more bold style
elems_only = datashader.spread(img)

nodes_only = datashader.dynspread(datashader.rasterize(trimesh.nodes).opts(opts.Image(cmap=['blue'])), shape='circle', max_px=6)

full_mesh = elems_only*nodes_only

elems_only.opts(alpha=0.2)*hv.Points(dfs, kdims=['x','y'], 
                                     vdims=['z','id','subloc','name']).opts(color='red', 
                                                                        size=10,
                                                                        tools=['hover'])
Out[8]:

Read and plot stations info¶

In [9]:
# param.nml has the time in the OPT section, but no parser for .nml files found

from datetime import datetime

import pandas as pd
import hvplot.pandas

def read_and_plot(file, station_file, reftime):
    df1 = station.read_staout(file, station_file, reftime)
    df1.index.name='Time' # workaround for hvplot bug
    return df1.hvplot()

reftime = datetime(2000,1,1)

plots = []
for index in range(1,10):
    fpath = f'../tests/data/m1_hello_schism/outputs/staout_{index}'
    station_file = '../tests/data/m1_hello_schism/station.in'
    vartype = schimpy.station.station_variables[index-1]
    try:
        plot = read_and_plot(fpath, station_file, reftime).opts(ylabel=f'{vartype}')
        plots.append(plot)
    except:
        pass
        #print(f'No data for index: {index}: {vartype}')

hv.Layout(plots).cols(1).opts(shared_axes=False)
Out[9]:

Select station points and display station out info¶

In [10]:
selectable_pts = hv.Points(dfs, kdims=['x','y'], 
                                     vdims=['z','id','subloc','name']).opts(color='red', 
                                                                        size=10,
                                                                        tools=['hover','tap'])

select_station = streams.Selection1D(source=selectable_pts)


def read_station_data(station_id, vartype, file_path_prefix, reftime):
    station_index = schimpy.station.station_variables.index(vartype)
    print(station_index)
    file = f'{file_path_prefix}{station_index+1}'
    df = station.read_staout(file, station_file, reftime)
    df.index.name='Time' # workaround for hvplot bug
    selected_cols = df.columns[df.columns.str.contains(station_id)]
    return df[list(selected_cols)]

def show_station_data(station_id, vartype, file_path_prefix, reftime):
    df = read_station_data(station_id, vartype, file_path_prefix, reftime)
    return df.hvplot()
 

#show_station_data('ocean', 'elev', '../tests/data/m1_hello_schism/outputs/staout_', reftime)

def show_station_data_for_index(index, vartype):
    if not index:
        return hv.Div('Select a point from the map')
    sdfs = dfs.iloc[index]
    ids = sdfs['id'].unique()
    return hv.Overlay([show_station_data(id, 
                                 vartype, 
                                 '../tests/data/m1_hello_schism/outputs/staout_',
                                 reftime) for id in ids])
In [11]:
vartype = pn.widgets.Select(options=station.station_variables)

pn.Column(elems_only.opts(alpha=0.2)*selectable_pts, 
          vartype,
           pn.Row(pn.bind(show_station_data_for_index, 
                          index = select_station.param.index, 
                          vartype = vartype)))
5
Out[11]:

Closest element ID on tap¶

In [12]:
import numpy as np

points = hv.Points([])
taps = hv.streams.Tap(source=points, x=np.nan, y=np.nan)

def location(x, y):
    if np.isnan(x) or np.isnan(y):
        element_id = np.nan
    else:
        element_id = grid.find_closest_elems([x,y])
    return pn.pane.Str(f'Click at {x:.2f}, {y:.2f}\nclosest element {element_id}')

pn.Column(points*elems_only, pn.bind(location, x=taps.param.x, y= taps.param.y))
Out[12]:
In [ ]: